__kernel void project(__global const float* p, __global float* result, const int rows) {
    const float eps = 1e-15;
    int i = get_global_id(0); // Each work-item handles one row of the matrix

    if (i < rows) {
        // Calculate matrix indices based on column-major order
        float x = p[i];
        float y = p[i + rows];
        float z = p[i + 2*rows];

        float r = sqrt(x * x + y * y);
        float s = atan2(r, z) / fmax(r, eps);

        result[i] = x * s;
        result[i + rows] = y * s;
    }
}